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Flowing blood displays a phenomenon called margination, in which leukocytes and platelets are 
preferentially found near blood vessel walls, while erythrocytes are depleted from these regions. 
Here margination is investigated using direct hydrodynamic simulations of a binary suspension of 
stiff (s) and floppy (f) capsules, as well as a stochastic model that incorporates the key particle 
transport mechanisms in suspensions - wall-induced hydrodynamic migration and shear-induced 
pair collisions. The stochastic model allows the relative importance of these two mechanisms to be 
directly evaluated and thereby indicates that margination, at least in the dilute case, is largely due 
to the differential dynamics of homogeneous (e.g. s-s) and heterogeneous (s-f) collisions. 



Introduction — Blood is a multicomponent suspension 
consisting primarily of red blood cells (RBCs) along with 
trace amounts of other components, primarily leukocytes 
and platelets [1, 2]. Under physiological flow conditions, 
both the leukocytes and platelets segregate near the ves- 
sel walls [3] , a phenomenon known as margination, while 
the RBCs tend to be depleted in the near-wall region 
forming a so-called "cell- free layer" [1, 2]. Leukocytes 
are larger than RBCs and platelets smaller, but both are 
considerably stiffcr than RBCs [4]; this difference is be- 
lieved to play an important role in their margination. In 
addition, in illnesses such as malaria and sickle cell dis- 
ease, the RBCs themselves are known to become stiff, 
and these stiffened cells also marginatc [5] . Furthermore, 
delivery of drugs to tumors using particles injected into 
the bloodstream is a major goal of current cancer re- 
search [6]; margination may influence the distribution 
and thus the efficacy of these particles. Finally, there 
are many microfluidic applications where differences in 
the margination properties of various blood components 
are exploited to effect their separation [2, 5, 7]. 

Despite the importance of flow induced segregation and 
margination phenomena in particle mixtures like blood, 
a mechanistic understanding is elusive [2]. The goals 
of the present work are to establish the mechanisms of 
rigidity-based margination in confined flows and to illus- 
trate that these mechanisms are generic for multicompo- 
nent suspensions. We do this in two parts. First, we con- 
sider direct hydrodynamic simulations of a model prob- 
lem that isolates the effect of stiffness on margination: 
a binary suspension of fluid-filled clastic capsules sub- 
jected to simple shear flow in a planar slit [8]. The two 
components of the binary mixtures have unequal mem- 
brane rigidities - the component with the higher rigidity 
is termed stiff, while the component with the lower rigid- 
ity is termed floppy. Additionally, we employ an idealized 
master equation (ME) model of the suspension dynam- 
ics that incorporates the two key transport mechanisms 
in confined suspensions: (1) wall-induced particle migra- 
tion and (2) particle displacements in homogeneous (e.g. 



stiff-stiff) and heterogeneous (stiff-floppy) pair collisions. 
We introduce a novel hydrodynamic Monte Carlo (HMC) 
simulation technique to find steady state concentration 
distributions for this model. In contrast to direct nu- 
merical simulations, the ME/HMC approach allows the 
various aspects of the particle dynamics to be indepen- 
dently controlled, thereby enabling delineation of their 
role in the segregation behavior. 

Using these approaches, we are able to isolate the ef- 
fect of heterogeneous pair collisions and thus demonstrate 
their dominant role, at least in dilute systems, in the ob- 
served margination behavior. The approach and mech- 
anisms presented here are generic for multicomponent 
suspensions - the model only takes as inputs the migra- 
tion behavior of the various components of the suspension 
and their displacements upon homogeneous and hetero- 
geneous collisions. Therefore, it is extensible to other 
systems including mixtures of particles of unequal sizes 
and shapes, thus encompassing whole blood flow as well 
as flow of drug delivery particles in the bloodstream. 

Formulation — We consider a suspension of fluid-filled 
nco-Hookean capsules [8] subjected to simple shear flow 
with shear rate 7 between infinite parallel walls at y = 
and y — H (Fig. la). At rest, all capsules are spheres 
with the same radius a. The rigidity of a particle is char- 
acterized by its membrane shear modulus G, which is ex- 
pressed in terms of the non-dimensional capillary number 
Ca = fxja/G, where \i is the viscosity of the suspending 
fluid. The capsules in a binary mixture with a lower 
Ca = Ca s are termed stiff, while the capsules with a 
higher Ca = Caj are termed floppy; the number fraction 
of the floppy particles will be denoted Xf. The viscosity 
ratio A of the fluid inside and outside the particle is unity. 
The Reynolds number is taken to be negligible. 

In a dilute suspension of particles, particle interactions 
can be treated as a sequence of uncorrelated pair col- 
lisions [9, 10]. In addition, since the capsules are de- 
formable, they also have a wall-induced migration ve- 
locity v m {y) away from the wall [11]. The effect of the 
pair collisions and the wall-induced particle migration 



2 



^ ^ ^ — s. 

/io°o 0d ^ ^ 
^ o o i 




y 
If 




FIG. 1. (color online) (a) System geometry, (b) schematic 
of particle trajectories in a pair collision, (c) isolated particle 
migration velocity at confinement ratio 2a/H — 0.f97, and 
(d) cross-stream displacement A in various types of pair col- 
lisions as a function of the initial offset 5. Here 's' refers to 
the stiffer particle (Ca a = 0.2), while 'f refers to the floppier 
particle (Ca/ = 0.5). 



can be consistently described by a kinetic master equa- 
tion (cf. [9]). For simplicity, we focus here on a model of 
a monolayer in the x — y (flow-gradient) plane - exten- 
sion to three dimensions is straightforward [9] . The mean 
area number density of all the particles in the monolayer 
is denoted no, such that the areal fraction in the x — y 
plane at rest is <fi a — ira 2 no. The mean area number den- 
sity of each of the species a in the mixture is denoted n$ , 
while its distribution in the y direction is denoted n a (y). 
In this case, the master equation is: 



dn a (y) _ S(ifin a ) 



dt 



dy 



N s , y 

■El/ [" Q fe- A<: 



x n f3 {y - A af} - S) - n a {y) n?(y - 6)] j \S\ dS^j , (1) 

where S is the pre-collision pair offset in the y direction, 
A a P(S) is the cross-stream displacement of particle of 
type a after collision with another particle of type /?, 
while the sum is over all the species N s in the suspension; 
see (Fig. lb) for a schematic of a pair collision. The first 
term on the right hand side arises from the wall-induced 
migration, while the integral term represents the effect 
of pair collisions [9]. 

This equation is analogous to the Boltzmann equation 
for rarefied gases [12, 13]. The dynamic simulation Monte 
Carlo (DSMC) approach is a popular technique to obtain 
solutions of the Boltzmann equation [12, 13]. The current 
work is inspired by the DSMC approach, which has also 
found interest in recent works on colloidal suspensions of 
rigid spheres [9] . By analogy, we term the method in the 



present work as the hydrodynamic Monte Carlo (HMC) 
method. As in the case of DSMC [12, 14], the HMC ap- 
proach is appropriate in the dilute limit and requires the 
assumption of chaotic particle dynamics, an assumption 
that is valid for the particulate flows considered here [15]. 

In the HMC approach, the y-positions of N p particles 
are followed in time. We set N p = 100 here; simulations 
with larger N p gave indistinguishable results. Each parti- 
cle is assumed to represent an infinite number of particles 
at the same y position randomly distributed in the flow 
direction with an average spacing of L, where L is given 
by L = N p / {tiqH). A distinguishing feature of the HMC 
(or DSMC) approach is that the collisions between par- 
ticles are treated probabilistically [14]. In the present 
study, we neglect pair collisions with large initial offsets 
5 > 5 cut , as their effect on cross-stream displacement is 
weak; we take 8 cu t = 2.5a here. A timestep of the simu- 
lation involves choosing a pair of particles which satisfies 
the condition 8 < S cu t- The pair is subsequently selected 
or rejected for collision with a probability proportional 
to their relative velocity of approach j\S\ [9, 12, 16]. An 
important aspect of the simulation is the time interval 
between collisions At, because the wall-induced parti- 
cle migration occurs simultaneously with the collisions. 
In order to determine At at each time step, we assume 
that the number of particle collisions with time follows a 
Poisson process with a mean collision frequency v [y is 
estimated [17]), such that the time interval between col- 
lisions is distributed with probability P(At) = ve 
[16]. The distribution of time interval between collisions 
results from the variation of particle positions in the flow 
direction, which is not explicitly specified. Once the col- 
lision pair and time interval are set, we update positions 
of all particles A; as: y^{t + At) = yk{t) + v m At + A&, 
where A^ is non-zero only for the colliding particle pair 
If the collision is homogeneous, then A^ = A ss 
or A^ for stiff and floppy particles, respectively. If the 
collision is heterogeneous, then A& = A s ? for the stiff 
particle and Afc = A^ s for the floppy one. The procedure 
outlined above is repeated until a statistical steady state 
is obtained. Other details of the method can be found in 
Koura [16], whose approach is closely followed here. 

Pair collisions and wall induced migration. — The 
HMC method requires as inputs the cross-stream dis- 
placements A"* 3 in pair collisions and the wall-induced 
migration velocity v^. These were computed using an 
accelerated boundary integral method [8, 18]; also see 
[19]. Figure (lc) shows the isolated particle migration ve- 
locity as a function of y/H for Ca s = 0.2 and Ca/ = 0.5 
capsules at a confinement ratio 2a/ H = 0.197. Results 
for A a/3 (S) for the same two species (in an unbounded 
domain) are shown in Fig. (Id). For offsets S < a, 
/\ff < Furthermore, an important feature in this 

plot, first reported in Ref. [8], is that the displacement 
of the stiffer particle in heterogeneous collisions is higher 
than that of the floppy particle (A s ^ > A-^ s ), while the 
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FIG. 2. (color online) (a) and (b): Mean normalized distance 
of a species d = 2d/H from the centerline in HMC (H) and BI 
(B) methods in (a) (Ca s , Ca/, <j> v ) = (0.2, 0.5, 0.04) mixture as 
a function of Xf and in (b) (0.2,0.5,0.12) mixture, (c) and 
(d): Normalized number density profile n for (c) stiff and (d) 
floppy particles at Xf — 0.5 for the suspension in (a). In all 
cases 2a/ H = 0.197. 

cross-stream displacements in homogeneous collisions are 
between these two limits, i.e., A fs < A ff , A ss < A< This 
ordering will turn out to be crucial in determining the 
segregation behavior. 

The above particle migration velocity and pair collision 
results were determined in idealized systems, namely, by 
considering an isolated particle and an unconfmed sys- 
tem, respectively. Our interest is in a confined suspen- 
sion, so we expect that corrections will be necessary to 
both the migration and the pair collision results. For 
the migration velocity v m , we note the wall-induced mi- 
gration is due to the disturbance velocity created by the 
particle's image, and is a far-field effect [11]. As a re- 
sult of its far-field nature, in a suspension of particles, 
it can be expected that the wall-induced migration of a 
particle will result not only due to interaction with its 
own images, but also due to images of other particles - 
this can be expected to introduce an averaging effect in a 
suspension of particle mixtures. To model this, we mod- 
ify the migration velocities as per the following equation: 

4 s + (1 - J2i s x P v m> where v m s is the mi- 

gration velocity of the species a in the suspension and £ 

is an adjustable parameter. By taking £ = 1, we recover 

the isolated particle migration velocity for each species, 

while for £ = 0, each of the species in the suspension has 

the same migration velocity. 

We adopt a similarly simple model to account for the 
confinement effects on cross-stream displacement A in 
pair collisions. For a spherical particle with center at 
y = a, its A in a pair collision is expected to be zero as it 
will be touching the wall, while its A will approach the 
unconfmed result at large particle-wall separations. To 
account for this effect, we multiply the A of the particle 
in an unconfmed system by a factor t] = 1 — e ~^ dm ~ a ^ a , 
where d w is the distance of the particle from the nearest 
wall, assumed to satisfy d w > a. Despite this correction 
factor, particles might still occasionally overlap with the 
walls at high volume fractions, though that did not occur 
in the regime investigated here. 

Validation. — We next validate the HMC method by 



comparing it with the results of detailed boundary in- 
tegral (BI) simulations [8, 18]. The suspensions in the 
BI simulations are fully three dimensional, but otherwise 
are similar to those in the HMC method. The particle 
volume fraction in the BI simulations will be denoted <f> v . 
The HMC model has two adjustable parameters, n (or 
4> a ) and £, which can be tuned to obtain good agreement 
with BI simulations for a given suspension. In the present 
work, we seek the agreement of the mean normalized dis- 
tance of a species from the centerline d = 2d/H (see Fig. 
la) between the two methods. 

The parameters no and £ in the HMC method are ex- 
pected to depend on <fi v and 2a/ H, while their depen- 
dence on Ca s , Ca/ and Xf is expected to be weak. To 
demonstrate this, we focus on suspensions with <fi v = 0.04 
and 2a/H = 0.197 held fixed. We then consider the BI 
results for a pure suspension with Ca = 0.2, and tune the 
value of no in the HMC method to obtain a good match 
in d between the two methods; this yields no = 0.026a~ 2 
(or (p a = 0.082). We next consider BI results for a bi- 
nary mixture with Ca s = 0.2, Ca/ = 0.5, and Xf = 0.5 
and tune the value of £ in the HMC method to obtain a 
good agreement in d of both the species between the two 
methods; this yields £ = 0.23. Having set the value of 
no and £, the HMC method can then be used to predict 
results for other suspensions at the same <f> v and 2a/ H. 
To establish this, we consider the same binary mixture 
as above (Ca s = 0.2, Ca/ = 0.5) and predict d for both 
the species for a range of Xf and compare them with the 
corresponding BI results (Fig. 2a). Excellent agreement 
can be observed at all values of Xf. Similar close agree- 
ment was also observed for different sets of (Ca s ,Ca/), 
namely (0.1, 0.5) and (0.3, 0.4) with no adjustment of no 
and £ (see [19]). Lastly, we consider a suspension with 
(Ca S! Ca/) of (0.2, 0.5) at cj> v = 0.12 and 2a/ H = 0.197, 
and we tune the values of no and £ as discussed above, 
which yield n = 0.093a" 2 (or 4> a = 0.292) and £ = 0.63. 
Subsequently, we predict d at various Xf in this system 
(Fig. 2b). Excellent agreement is observed even at this 
higher volume fraction. 

Besides the averaged measure d, we also compared 
the particle number density distribution in the wall nor- 
mal direction in the HMC and BI simulations. Re- 
sults for the normalized number density distribution 
n a (y) = n a (y)/no are shown in Figs. (2c) and (2d) for 
the Ca s = 0.2 and Ca/ = 0.5 mixture at cf> v = 0.04, 
2a/H = 0.197, and Xf = 0.5; a more complete set of 
plots is provided in [19]. Very good agreement of n with 
the BI results is observed in all cases. The agreement 
is remarkably good in the region around the centerline, 
though the peak near the wall is usually smeared in the 
HMC results in comparison to the BI results. Nonethe- 
less, given the broad agreement of d as well n in HMC 
and BI simulations in various suspensions, it is appar- 
ent that the HMC model captures the key aspects of the 
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FIG. 3. (color online) Number density h for stiff (top row) and 
floppy (bottom row) particles in simulations at Xf = 0.01 and 
Xf — 0.99, from: (a) and (d), full model (case A); (b) and (e), 
only difference between the two species being heterogeneous 
collisions (case B); (c) and (f), only difference between the 
two species being migration velocities (case C). These plots 
are for (Ca s , Ca/, <f> v ) = (0.2, 0.5, 0.04). 



particle distributions in these suspensions. 

Mechanisms of flow induced segregation. — A key ben- 
efit of the HM C approach is that it allows independent 
investigation of the effect of various ingredients of the 
particle dynamics on the margination behavior. We fo- 
cus here on suspensions with Ca s = 0.2, Ca/ = 0.5, 
<j) v = 0.04, and 2a/ H = 0.197. As discussed above, the 
HMC results for this system were generated by setting 
no = 0.026a -2 and £ = 0.23. We first show the plots for 
n from the full HMC model for both the stiff and floppy 
particles at Xf = 0.01 (dilute in floppy) and Xf = 0.99 
(dilute in stiff) in Figs. (3a) and (3d), respectively. It is 
clear from these plots that the stiff particles accumulate 
in the particle layer formed nearest to the wall as they be- 
come dilute in the suspension, i.e., they marginate with 
increasing Xf. In contrast, the floppy particles do the 
opposite as they become dilute (Xf decreases), accumu- 
lating near the centerline and thus "antimarginating" . 
These trends agree with the detailed boundary integral 
results of Ref. [8]; also see [19]. 

To disentangle the effects of wall induced migration 
and pair collisions, we now consider a number of control 
cases. First, we investigate the impact of heterogeneous 
collisions, by (i) setting the particle migration velocities 
of both the species to the simple average migration ve- 
locity of these two species and (ii) setting A"" and to 
the average value for the two species. Therefore, the only 
difference between these two species is their behavior in 
heterogeneous collisions: A"^ > A^ s . Plots of h(y) for 
these simulations are shown in Figs. (3b) and (3e). The 
difference between A s f and A? 3 is sufficient to lead to a 
segregation between the two species. In fact, as will be 
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FIG. 4. (color online) The difference Ad between d of a 
species in the mixture and the d of that species in its pure 
suspension for (a) stiff particles and (b) floppy particles in 
(Ca s , Ca/, cj> v )=(0.2, 0.5, 0.04) suspension. The cases A-C are 
described in the text. 



quantified shortly, most of the segregation results from 
heterogenous collisions. 

We next examine the effect of differences in migration 
velocity on the segregation behavior by setting the cross- 
stream displacement in all types of collisions for both the 
species to the simple average of the four curves on Fig. 
(Id), yielding A sf = A ss = A? 3 = A#. Plots for h(y) in 
this case are shown in Figs. (3c) and (3f). Here too some 
segregation is observed, though the degree of segregation 
is considerably smaller than in the full model. 

The degree of segregation between the two species is 
more quantitatively characterized by computing the dif- 
ference in d of each of the species from the corresponding 
pure species result; this is denoted by Ad. The plots for 
Ad in various cases described above are shown in Figs 
(4a) and (4b) for the stiff and floppy particles, respec- 
tively. For the present parameter set, the degree of segre- 
gation from the full model (case A) and the model where 
only heterogeneous collisions are distinct (case B) are al- 
most identical, while that resulting from differences in the 
migration velocity (case C) is much weaker. In a recent 
direct simulation study [8] , heterogeneous collisions were 
conjectured to play an important role in margination. In 
the simulations, however, it is not possible to deconvolve 
the effects of migration and collisions: i.e. one cannot 
perform control simulations. The ME/HMC approach 
does not suffer from this limitation. 

Finally, we have investigated the effect of volume frac- 
tion on particle dynamics and margination. At <p v = 
0.12, the effects of heterogeneous collisions and differen- 
tial migration velocities (cases B and C) on margination 
become comparable (see [19]). Very similar results arise 
at (j) v = 0.2 [19]. In the microcirculation, the volume 
fraction of cells is between <p v = 0.1 — 0.25 and it appears 
that in this regime, both the migration velocity and the 
heterogeneous collisions are playing an important role in 
the segregation between the various species. 

Conclusions. — To gain a mechanistic understanding 
of margination in blood and other multicomponent sus- 
pensions, we have introduced an idealized description 
of the particle flow dynamics that incorporates the two 
key transport mechanisms in confined suspensions: wall- 



5 



induced migration and hydro-dynamic pair collisions. The 
results for this model system clarify the important and 
previously underappreciated role played by heteroge- 
neous collisions in the observed segregation behavior. Be- 
cause differential behavior in heterogeneous collisions is 
generic for particles with contrasting shape or size as well 
as flexibility, the insights presented here are also applica- 
ble for other multicomponent suspensions. In particular, 
they have important implications in the design of drug 
delivery particles for optimal vascular wall targeting or 
for separating trace components of blood in microfluidic 
devices. 
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